global  uint128_div_v1_asm

BITS 64

SECTION .text  align=16

; ********** uint128_div_v1_asm function **********
; C/C++ prototype:
; extern "C" void uint128_t uint128_div_v1_asm (uint128_t n, uint128_t d);
; Input:   n, d,
; Output:  return n/d

; input parameter
; n:  RSI:RDI,  RDI: low64(n), RSI: high64(n)
; d:  RCX:rdx,  RDX: low64(d), RCX: high64(d)
; ***********************************************************

uint128_div_v1_asm :

%define N_HI rsi
%define N_LO rdi
%define D_HI r10
%define D_LO r11
%define T_HI r8    ; tmp_hi
%define T_LO r9    ; tmp_low
%define CNT  rcx ; counter

    mov D_LO, rdx           ; save n to D_HI:D_LO
    mov D_HI, rcx

    mov rax, N_LO
    mov rdx, N_HI
    sub rax, D_LO
    sbb rdx, D_HI           ; rdx:rax = n-d

    jc  return_ZERO         ; if n<d goto return_ZERO

    or D_HI, D_HI
    jz D_HI_IS_0            ; n is "xx" form, do xx/0y or 0x/0y
    jmp   L_XX_DIV_YY_IMPL  ; n="xx" form and d is 'yy" form, do xx/yy

D_HI_IS_0:                  ;  d<2^64
    or D_LO,D_LO
    jz DIV_OVERFLOW

    or  N_HI, N_HI
    jz  L_0X_DIV_0Y        ; if n is "0x" form, then do 0x/0y
    jmp L_XX_DIV_0Y_IMPL   ; n="xx" form and d is '0y" form, do xx/y

L_0X_DIV_0Y:               ; n="0x" form, n<2^64 and d<2^64
   jmp  L_0X_DIV_0Y_IMPL

;-----------------------------------------
; *** case 1:  n < d , return 0 **
align 16
return_ZERO:
    xor rdx,rdx
    xor rax,rax
    ret

;-----------------------------------------
; *** case 2:  n/0
align 16
DIV_OVERFLOW:   ; do 1/0
    mov rax, 1
    div D_LO
    ret

;-----------------------------------------
; *********  case 3:  do xx / yy *********
align 16
L_XX_DIV_YY_IMPL:
    bsr	CNT, N_HI   ; CNT=bit count of N_HI
    bsr	rax, D_HI   ; CNT=bit count of D_HI
    sub CNT, rax    ; shift = log2(n)-log2(d)

    shld D_HI, D_LO, cl
    shl	 D_LO, cl   ; D_HI:D_LO = (d<<shift)
    add CNT, 1      ; shift++

%undef R_MASK
%undef R_Q
%define R_MASK rdx
%define R_Q    rax

    xor R_Q,R_Q             ; q=0
    ; because n>=d, so CNT>=1, always enter loop

 align 16
loop1_start:
    xor R_MASK, R_MASK      ; mask=0
    sub N_LO, D_LO
    sbb N_HI, D_HI          ;  n-=d
    sbb R_MASK, 0           ;  if (n<0),then mask=-1, else mask=0

    lea R_Q, [R_Q*2+R_MASK+1]  ; final q=q*2+(mask+1)

    mov T_LO, D_LO          ; tmp_low =  low64(d)
    mov T_HI, D_HI          ; tmp_high = hi64(d)
    and T_LO, R_MASK        ; if mask==-1, then T_LO=D_LO, else T_LO=0
    and T_HI, R_MASK        ; if mask==-1, then T_HI=D_HI, elst T_HI=0

    add N_LO,T_LO
    adc N_HI,T_HI           ;  restor N

    shrd D_LO, D_HI, 1
	shr	 D_HI, 1            ; d = d/2

    sub  CNT,1              ; shift--
    jnz loop1_start         ; if shift>0, loop again

    xor rdx, rdx            ; if shift==0, return 0:R_Q, R_Q is rax
    ret

;-----------------------------------------
; ******  case 4: do X/Y ****
L_0X_DIV_0Y_IMPL:
    bsr	CNT, N_LO
    bsr	rax, D_LO
    sub CNT, rax        ; shift = log2(n)-log2(d)
    shl D_LO, cl        ; d= (d<<shift)
    add CNT,1           ; shift ++
%undef R_MASK
%undef R_Q
%define R_MASK rdx
%define R_Q    rax

    xor R_Q,R_Q             ; q=0
    ; because n>=d, so CNT>=1, always enter loop

    align 16
loop2_start:
    xor R_MASK, R_MASK      ; mask=0
    sub N_LO, D_LO          ; n-=d
    sbb R_MASK, 0           ; if (n<0),then mask=-1, else mask=0

    lea R_Q, [R_Q*2+R_MASK+1]  ; final: q=q*2+(mask+1)

    and R_MASK, D_LO        ; if mask==-1, then mask=D_LO, else mask=0
    add N_LO,   R_MASK      ; restor N

    shr  D_LO,  1           ; d = d/2

    sub  CNT,1              ; shift --
    jnz loop2_start         ; if shift>0, loop again

    xor rdx, rdx            ; R_Q is rax
    ret

 ;-----------------------------------------
 ; *********        case 5: do xx/Y     ******
align 16
L_XX_DIV_0Y_IMPL:
    bsr	CNT, N_HI   ; count= log2(n)-64
    bsr	rax, D_LO   ; rax= log2(d)
    add CNT, 64     ; count=  log2(n)
    sub CNT, rax    ; shift = log2(n)-log2(d)

    xor rax, rax        ; rax=0
    shld D_HI, D_LO, cl
    sal	 D_LO, cl       ; d= (d<<shift)

    test    CNT,  64
    cmovne	D_HI,  D_LO
    cmovne	D_LO,  rax

    add CNT, 1

%undef T_HI
%undef R_MASK
%undef R_Q
%define R_MASK r8
%define Q_HI rdx
%define Q_LO rax

    xor Q_LO,Q_LO
    xor Q_HI,Q_HI           ; q=0
    ; because n>=d, so CNT>=1, always enter loop

    align 16
loop3_start:
    xor R_MASK, R_MASK      ; mask=0

    sub N_LO, D_LO
    sbb N_HI, D_HI          ; n-=d
    sbb R_MASK, 0           ; if (n<0),then mask=-1, else mask=0

    add Q_LO, Q_LO
    adc Q_HI, Q_HI          ; q *= 2

    lea Q_LO, [Q_LO+R_MASK+1]    ; final: q = q*2 + mask+1

    mov T_LO, D_LO          ; tmp_low=  low64(d)
    and T_LO, R_MASK        ; if mask==-1, then T_LO=D_LO,   else T_LO=0
    and R_MASK, D_HI        ; if mask==-1, then mask =D_HI,  else mask=0

    add N_LO,T_LO
    adc N_HI,R_MASK         ; restor N

    shrd D_LO, D_HI, 1
	shr	 D_HI, 1            ; d = d/2

    sub  CNT,1              ; shift--
    jnz loop3_start         ; if shift>0, loop again

    ret                     ; Q_LO is rax, and Q_HI is rdx, do nothing
